Importing GPS data into R

require(raster);require(sp);require(rgdal); require(ggmap)
## Loading required package: raster
## Loading required package: sp
## Loading required package: rgdal
## rgdal: version: 1.2-18, (SVN revision 718)
##  Geospatial Data Abstraction Library extensions to R successfully loaded
##  Loaded GDAL runtime: GDAL 1.11.3, released 2015/09/16
##  Path to GDAL shared files: /usr/share/gdal/1.11
##  GDAL binary built with GEOS: TRUE 
##  Loaded PROJ.4 runtime: Rel. 4.9.2, 08 September 2015, [PJ_VERSION: 492]
##  Path to PROJ.4 shared files: (autodetected)
##  Linking to sp version: 1.2-7
## Loading required package: ggmap
## Loading required package: ggplot2

#  Set your working directory
setwd('/home/pjg/GIS')
# Create a directory folder and move into it
dir.create('AdvancedGIS_R')
## Warning in dir.create("AdvancedGIS_R"): 'AdvancedGIS_R' already exists
setwd('AdvancedGIS_R')

If you are using data from a Garmin GPS

##############################################################
###############  Importing Garmin GPS data  ##################
##############################################################
###### Create a character vector of the entire pathway for the gpx object from your Garmin GPS
gpfile<-"/home/pjg/GIS/Classes/spatial_bioinformatics/R_GIS/Current.gpx"
###### Read in the GPX file using the function from rdgal
GPSfile<-readOGR(dsn = gpfile, layer = "track_points")
## OGR data source with driver: GPX 
## Source: "/home/pjg/GIS/Classes/spatial_bioinformatics/R_GIS/Current.gpx", layer: "track_points"
## with 3 features
## It has 26 fields
##### Extract all of the data and coordinates
wpsDF<-GPSfile@data
wpsCoords<-GPSfile@coords
##### Pull out only the necessary data:longitude, latitude, time, elevation
##### Garmin does not supply accuracy at each waypoint
wps<-cbind(wpsCoords[,1:2], wpsDF[,c(5,4)])
colnames(wps)<-c("lon", "lat", "time", "elev")        
head(wps)
##         lon      lat                   time  elev
## 1 -73.94220 40.85085 2018/05/06 20:09:54+00 81.85
## 2 -73.94218 40.85087 2018/05/06 20:10:28+00 81.85
## 3 -73.94220 40.85085 2018/05/06 20:10:57+00 81.37

If you are using GPS data from Cellphone

###############################################################
############  Importing cellphone data  #######################
###############################################################
# Read in .csv file of waypoints
wps<-read.csv('/home/pjg/GIS/Classes/spatial_bioinformatics/R_GIS/20180413.csv')
# Take a look at the first few rows of the data.frame
head(wps)
##                       time      lat       lon elevation accuracy bearing
## 1 2018-04-13T15:01:31.966Z 41.22050 -74.29108        NA   20.274      NA
## 2 2018-04-13T15:03:20.109Z 41.22051 -74.29101        NA   20.249      NA
## 3 2018-04-13T15:03:20.110Z 41.22051 -74.29101        NA   20.249      NA
## 4 2018-04-13T15:04:47.000Z 41.22035 -74.29116       263   37.000      NA
## 5 2018-04-13T15:05:08.000Z 41.22019 -74.29070       254   39.000      NA
## 6 2018-04-13T15:05:36.119Z 41.22050 -74.29095        NA   20.460      NA
##   speed satellites provider hdop vdop pdop geoidheight ageofdgpsdata
## 1    NA          0  network   NA   NA   NA          NA            NA
## 2    NA          0  network   NA   NA   NA          NA            NA
## 3    NA          0  network   NA   NA   NA          NA            NA
## 4    NA         11      gps  1.2  0.9  1.5         -37            NA
## 5     0         11      gps  1.2  0.9  1.5         -37            NA
## 6    NA          0  network   NA   NA   NA          NA            NA
##   dgpsid activity battery annotation
## 1     NA       NA      90         NA
## 2     NA       NA      90         NA
## 3     NA       NA      90         NA
## 4     NA       NA      90         NA
## 5     NA       NA      90         NA
## 6     NA       NA      89         NA
# We are only really interested in the first 5 columns of the waypoint file, but let's rearrainge them a bit
wps<-wps[c(3,2,1,4,5)]
# Now lon and lat are the first 2 columns of the data.frame
head(wps)
##         lon      lat                     time elevation accuracy
## 1 -74.29108 41.22050 2018-04-13T15:01:31.966Z        NA   20.274
## 2 -74.29101 41.22051 2018-04-13T15:03:20.109Z        NA   20.249
## 3 -74.29101 41.22051 2018-04-13T15:03:20.110Z        NA   20.249
## 4 -74.29116 41.22035 2018-04-13T15:04:47.000Z       263   37.000
## 5 -74.29070 41.22019 2018-04-13T15:05:08.000Z       254   39.000
## 6 -74.29095 41.22050 2018-04-13T15:05:36.119Z        NA   20.460

Visualizing waypoints on satellite imagery

###############################################################
############  Visualizing and creating features  ##############
###############################################################
##  Visualizing your waypoints with satellite data
# Find the bounding box of your waypoints
loc = cbind(c(min(wps$lon), max(wps$lon)), c(min(wps$lat), max(wps$lat)))
# Name the columns
loc = as.data.frame(loc)
colnames(loc) = c('lon', 'lat')
# Get the images from Google Earth Engine
manbox <- make_bbox(lon = loc$lon, lat = loc$lat, f = .1)
manmap <- get_map(location = manbox, maptype = "satellite", source = "google", zoom =18)
# Plot
ggmap(manmap) + 
  geom_point(data = loc, 
             color = "red",
             size =1) # for better resolution, change size to smaller

Creating Shapefile Features from localities

##  Transform the waypoints and save as a shapefile
wpsShp<-SpatialPointsDataFrame(coords = wps[,1:2], data = wps)
#writeOGR(obj = wpsShp, dsn = paste(getwd()), layer = "Allwaypoints", driver = "ESRI Shapefile")

##  Separate the data into the appropriate features from which waypoints were taken outside

The first step is to take another look at the data and find the rownumber of the corresponding feature

print(wps)
##          lon      lat                     time elevation accuracy
## 1  -74.29108 41.22050 2018-04-13T15:01:31.966Z        NA   20.274
## 2  -74.29101 41.22051 2018-04-13T15:03:20.109Z        NA   20.249
## 3  -74.29101 41.22051 2018-04-13T15:03:20.110Z        NA   20.249
## 4  -74.29116 41.22035 2018-04-13T15:04:47.000Z       263   37.000
## 5  -74.29070 41.22019 2018-04-13T15:05:08.000Z       254   39.000
## 6  -74.29095 41.22050 2018-04-13T15:05:36.119Z        NA   20.460
## 7  -74.29109 41.22045 2018-04-13T15:06:36.470Z        NA   23.969
## 8  -74.29080 41.22039 2018-04-13T15:07:19.000Z       216   38.000
## 9  -74.29077 41.22048 2018-04-13T15:06:54.837Z        NA   21.048
## 10 -74.29083 41.22043 2018-04-13T15:07:36.000Z       204   37.000
## 11 -74.29064 41.22059 2018-04-13T19:06:07.000Z       180   21.000
## 12 -74.29064 41.22055 2018-04-13T19:07:38.000Z       186   40.000
## 13 -74.29066 41.22055 2018-04-13T19:08:47.000Z       188   34.000
## 14 -74.29067 41.22056 2018-04-13T19:10:09.000Z       186   31.000
## 15 -74.29068 41.22057 2018-04-13T19:10:13.000Z       185   16.000
## 16 -74.29068 41.22057 2018-04-13T19:10:16.000Z       185   18.000
## 17 -74.29120 41.22047 2018-04-13T19:09:43.811Z        NA   19.906
## 18 -74.29120 41.22047 2018-04-13T19:09:43.811Z        NA   19.906
## 19 -74.29069 41.22057 2018-04-13T19:10:25.000Z       185   23.000
## 20 -74.29070 41.22057 2018-04-13T19:10:28.000Z       184   16.000
## 21 -74.29071 41.22058 2018-04-13T19:10:30.000Z       183   15.000
## 22 -74.29115 41.22049 2018-04-13T19:10:07.017Z        NA   22.561
## 23 -74.29082 41.22050 2018-04-13T19:20:28.000Z       173   32.000
## 24 -74.29077 41.22047 2018-04-13T19:20:11.969Z        NA   21.443
## 25 -74.29078 41.22048 2018-04-13T19:20:26.408Z        NA   22.295
# Enter the row number of the point feature below
PointFeat<- 1
# Enter the row numbers of the line features below, separated with a colon
LineFeat<- 2:10
# As above, enter the row numbers of the polygon features below, separated with a colon  
PolyFeat<-  11:22 
# Now, use these row numbers to create individual 'spatial' objects for each of the features
# Point features
PointShape <- SpatialPointsDataFrame(wps[,1:2][PointFeat,], data = wps[PointFeat,])
# Line features
PointLine <- SpatialLinesDataFrame(sl = SpatialLines(list(Lines(list(Line(coords = wps[,1:2][LineFeat,])), ID='a'))), match.ID = F, data = as.data.frame("Line"))
# Polygon features
PointPoly <- SpatialPolygonsDataFrame(Sr = SpatialPolygons(list(Polygons(list(Polygon(wps[PolyFeat,][,1:2])),"Poly"))), data = data.frame("Poly"), match.ID = F)

Saving as ESRI shapefiles

writeOGR(obj = PointShape, dsn = paste(getwd()), layer = "PointFeature", driver = "ESRI Shapefile")
writeOGR(obj = PointLine, dsn = paste(getwd()), layer = "LineFeature", driver = "ESRI Shapefile")
writeOGR(obj = PointPoly, dsn = paste(getwd()), layer = "PolyFeature", driver = "ESRI Shapefile")

Test for accuracy

PointShape
ggmap(manmap) + 
  geom_point(data = data.frame(PointShape@coords), 
             color = "red",
             size =1) # for better resolution, change size to smaller

# LineShape
ggmap(manmap) + 
  geom_polygon(data = PointLine, 
               aes(x = long, y = lat),
               color = "red",
               size =1) # for better resolution, change size to smaller

# PolyShape
ggmap(manmap) + 
  geom_polygon(data = PointPoly, 
               aes(x = long, y = lat),
               color = "red",
               size =1) # for better resolution, change size to smaller
## Regions defined for each Polygons